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INTRODUCTION 

Apollo landing videos shot from inside the right LEM window, provide a 
quantitative measure of the characteristics and dynamics of the ejecta spray of 
lunar regolith particles beneath the Lander during the final 1 0 [m] or so of descent. 
Photogrammetry analysis gives an estimate of the thickness of the dust layer and 
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Figure 1. 2D CFD boundary and grid point definition. 



angle of trajectory (Immer, 2008). In addition, Apollo landing video analysis 
divulges valuable information on the regolith ejecta interactions with lunar surface 
topography. For example, dense dust streaks are seen to originate at the outer rims 
of craters within a critical radius of the Lander during descent. 

The primary intent of this work was to develop a mathematical model and 
software implementation for the trajectory simulation of lunar dust particles acted 
on by gas jets originating from the nozzle of a lunar Lander, where the particle 
sizes typically range from 10 pm to 500 pm. The high temperature, supersonic jet 
of gas that is exhausted from a rocket engine can propel dust, soil, gravel, as well 
as small rocks to high velocities. The lunar vacuum allows ejected particles to 
travel great distances unimpeded, and in the case of smaller particles, escape 
velocities may be reached. The particle size distributions and kinetic energies of 
ejected particles can lead to damage to the landing spacecraft or to other hardware 
that has previously been deployed in the vicinity. Thus the primary motivation 
behind this work is to seek a better understanding for the purpose of modeling and 
predicting the behavior of regolith dust particle trajectories during powered rocket 
descent and ascent. 



Figure 2. 3D DSMC simulation geometry. 


COMPUTATIONAL STRATEGY 


The trajectory model described in this paper requires output from a two- or three- 
dimensional fluid dynamics simulation of the rocket nozzle exhaust. Traditionally, 
this is achieved using a commercially available computational fluid dynamics 
(CFD) software package. Alternatively, direct simulation Monte Carlo (DSMC) 
software may be employed. DSMC simulators model fluid flows using a large 
number of discrete molecules in a Monte Carlo fashion to solve the Boltzmann gas 
equation. The worked described in this paper utilizes both CFD and DSMC 
simulations of an Apollo-like lunar Lander. In the case of the CFD simulations, 
2D symmetry was exploited based on cylindrical symmetry where the problem is 
independent of azimuth angle, reducing complexity and computation time (see 
Figure 1). Full 3D geometries were processed by DSMC software, as depicted in 
Figure 2, where no symmetries were assumed. 

The CFD/DSMC output provides estimates of gas density p(r) , gas velocity u(r), 
and gas temperature r(r) for a gridded volume described by r k vector points in the 

bounded domain. These values are interpolated from the CFD/DSMC grid by 
finding four nearest grid neighbors in a volume around the Ath trajectory point and 
applying an TV-dimensional interpolation algorithm. Inputs which are initial 
conditions of the trajectory calculation include the particle diameter D and initial 
position of the particle r 0 = (xo,yo, ^o) where the vertical direction x is positive up 
and equal to zero at the surface. Typically, the particle starting position above the 
surface is set by the user as x 0 . 

MATHEMATICAL MODEL AND TRAJECTORY ALGORITHM 

For a particle of diameter D and mass m, the trajectory is due to three external 
forces: lunar gravity, jet gas drag, and lift caused by the gas flow. The sum of 
external forces is equal to the acceleration experienced by the particle, which can 



be estimated by a Taylor series expansion about time point k, resulting in a set of 
difference equations for position and velocity: 
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g L is lunar gravity and p L is the lunar soil particle density (3100 [kg m' 3 ]);p(r) 
is the gas density and u(r) is the gas velocity from the CFD/DSMC files (*.g and 
*.q files). These values are computed by finding nearest neighbors in a volume 
around the point of interest and applying the spatial interpolation algorithm of 
Shepard (1968). 


a ufl is the lift acceleration due to the vertical gradient of the horizontal component 
of gas flow (Shao, 2000): 
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where C L is the coefficient of lift. 

The direction of the lunar gravity unit vector is given by: 
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where the lunar gravity constant g L is defined at the surface of a sphere of radius 
R l (moon radius). 

The coefficient of drag, C D is a function of the Reynolds number, R : 


R k = 
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where T k is the static temperature of the gas at the &th position of the particle,. The 
empirical parameters in the denominator are: A = 1.71575x1 O' 7 and /?= 0.78. 

The coefficient of drag is computed from the following empirical formula: 
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The initial conditions are: 


v 0 =0 
r o= R o 
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SIMULATION EXPERIMENTS 


Numerous gas jet simulations involving variations of rocket nozzle height above 
ground, nozzle angles, and number of craters at various distances and diameters, 
were run using DSMC software. Some of these variations are depicted in the 
schematic of Figure 2. The DSMC output then became the particle trajectory 
model (PTM) input, where PTM is a FORTRAN code implementation of 
Equations (1) through (6). These results were reported in detail by Metzger 
(2007). 

For the purpose of this present work, only CFD output is used as inputs to PTM. 
The four CFD test cases used in the trajectory simulation experiments are 
summarized in Tables 1. 


Table 1. CFD Cases used as Inputes to PTM. 



Distance from Nozzle Center Line to 
Crater Center, R c 

Height of Nozzle 
above Surface, h 

5 Rn 

15 Rn 

2.5 R n 

Case C2 

— 

5 Rn 

Case C7 

Case C3 

10 R n 

Case Cl 

— 


Figures 3, 4, and 5 display the CFD Case 2 output array values needed for input to 
PTM, i.e., gas density, gas velocity, and gas temperature. 


Table 2. Example PTM input parameters. 
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Table 3. PTM time step versus particle diameter. 


D [pm] 

10 

25 

50 

100 



At [ps] 

200 

500 

1000 

2000 

5000 

10000 


Figure 6 shows the PTM output for six particle sizes given by Table 3, spanning 
the range of lunar dust particle sizes. The initial starting coordinate in Table 2 is 
used in all six runs and corresponds to a point 3 [mm] above the outer lip of the 
crater. The reason for choosing a point slightly above the ground will be 
discussed in the following section. Similar plots are generated for the other cases 
shown in Table 1 . 

In order to compare results from the four CFD cases studied above, the velocity 
vector corresponding to the trajectories of each of the six particle sizes is 
calculated from the trajectory data. The horizontal right boundary at 10 [m] is a 
convenient point to evaluate velocity. Figure 8 is a particle velocity magnitude 
(speed) plot of three CFD cases, corresponding to three rocket nozzle heights, h, 
above the ground compared with DSMC data from Metzger (2007). Figure 9 is 
the particle velocity angles relative to the ground for the three CFD cases. 

PTM plots for six particle sizes is shown in Figure 7, corresponding to Case C3. 

In Case C3, the crater is located at 15^v * 30 [ft]. Initializing the starting point of 
the particles according to Table 2, results in velocities with lower angles and 
speeds. This result is a consequence of the turbulence generated above the crater, 
slowing down the particle with a downward component of velocity. 
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Figure 3. Example gas density plot corresponding to CFD Case 
C2. Horizontal and vertical axes are marked in meters. 



Figure 4. Example gas velocity plot corresponding to CFD 
Case C2. Horizontal and vertical axes are marked in meters. 



Figure 5. Example temperature contour plot corresponding to 
Case C2. Horizontal and vertical axes are marked in meters. 



Figure 6. PTM output for CFD Case 2 and initial values from 
Table 2. Horizontal and vertical axes are marked in meters. 




Figure 7. PTM output for CFD Case 3 and initial values from 
Table 2. Horizontal and vertical axes are marked in meters. 




Figure 8. Particle speeds exiting the CFD boundary for various particle sizes and 
CFD Cases of Table 1. The solid black line is averaged data from Metzger (2007). 



Figure 9. Particle trajectory angles relative to ground corresponding to 
Figure 7. Again, the solid black line is averaged data from Metzger (2007). 


DISCUSSION AND CONCLUSIONS 


The particle velocities converge nicely for sizes in the range of 200 to 500 [pm], 
as shown by Figure 8. For smaller particle sizes, there is a trend that shows higher 
velocities for larger nozzle height above the ground for the three cases considered, 
i.e., h = 5, 10, and 20 [ft]. The velocity magnitudes in this region are significantly 
larger than those shown by Metzger (2007). This may be attributed to the particle 
starting point: in the present case, particles were stated from the crater outer rim, 
whereas in Metzger’s data, the starting point was typically not on a crater rim. 

The CFD trajectory angle results shown in Figure 9 agree well with the 
photogrammetry results of Immer (2008). However, these results are not in good 
agreement with Metzger (2007), except for a narrow region around particle sizes 
of 200 [pm]. The source of this disagreement again may be due to differences in 
DSMC versus CFD simulations. It should also be noted that the coefficient of lift 
used in the current work is set to C L = 20, whereas the previous work by Metzger, 
C L = 500. 

In the absence of a vertical component of drag, which is a reasonable assumption 
very near the ground (if the ground is smooth), the lift force must overcome the 
weight of the particle in order to lift the particle into the gas stream where 
horizontal drag takes over and dominates particle trajectory: 
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where g L is lunar gravity (1.622 [m s' 2 ]), p L is the lunar soil particle density 
(3100 [kg m' 3 ]), u k is the particle speed relative to gas speed at time step k, and 


dx 


is the gradient of u k in the vertical direction. For example, using CFD Case 3 


with the initial particle coordinates given by Table 2, u 0 * 1 10 [ m s' 1 ] for a 50 


[pm] diameter particle, -^-*10 [s' 1 ], and p 0 *2x10 5 [igra' 3 ]. Evaluating 

Equation (7) with these values with Q = 20, a Lifi ~ 4.3 [m s 2 ] > g L . Therefore, 
the particle experiences lift. 

As a counter example, using the same initial conditions as above but let the 
vertical initial coordinate start with the particle resting on the ground, 

du 

D / 2 = 25 [pm] where D is the particle diameter: u 0 « 1 0 [ms 1 ], — - » 1 0 [ s' 1 ] , 

dx 

and again p 0 * 2 x 1 0" 5 [Ag m' 3 ] . Evaluating Equation (7) with these values with 
Ci = 20, a uft « 0.4 [m s' 2 ] <g L . In this case there is no lift. Note that using a 
larger value of Q = 500 (Metzger, 2007), lift is achieved in this example. 

The forces experienced by the particle are primarily due to drag, i.e., the particle 
will tend to follow the gas jet stream. However, particles on a perfectly smooth 
ground experience very little drag because of the no slip boundary condition which 
results in zero gas velocity at the ground boundary. Under these conditions, a lift 
force can artificially be generated in the trajectory simulation by using a very large 
value of lift coefficient (Metzger, 2007). More work will need to be done to better 
understand the conditions of the gas at the boundary. A hint of what might really 
be going on near the surface is observed by the disruption of flow caused by 
craters. On a smaller scale, any variations of the lunar surface are likely leading to 
boundary conditions that account for dust particle lift, beyond the simple model of 
aerodynamic lift presented above. Future work might also attempt to model the 
collision of particles leading to something like a chain-reaction, causing particles 
to pop up into the gas jet streams from the surface. Once in the gas flow, they are 
easily dragged into a trajectory and ejected at the approximate 3 degree angle 
observed in the Apollo videos. 
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